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Thermal density functional theory (DFT) calculations often use the Mermin-Kohn-Sham (MKS) 
scheme, but employ ground-state approximations to the exchange-correlation (XC) free energy. In 
the simplest solvable non-trivial model, an asymmetric Hubbard dimer, we calculate the exact many- 
body energies, the exact Mermin-Kohn-Sham functionals for this system, and extract the exact XC 
free energy. For moderate temperatures and weak correlation, we find this approximation to be 
excellent. We extract various exact free energy correlation components and the exact adiabatic 
connection formula. 
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I. INTRODUCTION 

Recent decades have seen enormous advances in the use 
of DFT calculations [1] of warm dense matter, a highly 
energetic phase of matter that shares properties of solids 
and plasmas [2]. Materials under the extreme tempera¬ 
tures and pressures necessary to generate WDM can be 
found in astronomical bodies, within inertial confinement 
fusion capsules, and during explosions and shock physics 
experiments [3]. These calculations are used in the de¬ 
scription of planetary cores [4-6], for the development of 
experimental standards [7, 8], for prediction of material 
properties [9-11], and in tandem with experiments push¬ 
ing the boundaries of accessible conditions [12]. Because 
of this growing interest in WDM and thermal systems 
in general, we seek to better understand thermal DFT 
using exactly solvable models. 

In almost all thermal DFT calculations, a crucial ap¬ 
proximation is made: the exchange-correlation (XC) free 
energy in principle depends on the temperature [13, 14], 
but in practice is approximated by a standard ground- 
state approximation. Most calculations are for extended 
systems, and usually use a generalized gradient approx¬ 
imation, such as PBE[15]. These Mermin-Kohn-Sham 
(MKS) [16, 17] calculations predict several key proper¬ 
ties, such as the free energy and density for a given dis¬ 
tribution of the nuclei, and any properties that can be ex¬ 
tracted from these, such as equations of state of materials 
and Hugoniot shock curves [18]. If the exact temperature- 
dependent XC free energy were known, such properties 
would be exact [19]. In some cases, response properties 
are extracted from the thermal KS orbitals [9], which in¬ 
volves a further approximation. Although no one has 
shown that the lack of thermal XC corrections is a fa¬ 
tal flaw in a given calculation, the pervasive use of this 
uncontrolled approximation is an underlying concern [20] 
that warrants investigation. 

The crucial step that made zero-temperature DFT suf¬ 
ficiently accurate for chemical purposes was the intro¬ 


duction and testing of generalized gradient approxima¬ 
tions about 20 years ago [21, 22]. By careful compari¬ 
son with highly accurate benchmarks produced either by 
direct solution of the Schrodinger equation or from ex¬ 
periments with well-controlled errors, the general level 
of accuracy and reliability of such approximations was 
well documented [23-25]. With improved binding ener¬ 
gies came the ability to determine molecular geometries 
for complex systems. A similar transformation is occur¬ 
ring in materials science today [26]. 

But no such database or highly accurate results exist 
for thermal systems. It is hard to imagine experimen¬ 
tal measurements of energies with the required accuracy 
under the relevant conditions, but calculations should be 
possible. Various Monte Carlo methods have been devel¬ 
oped to study WDM in extended systems [27-33]. There 
have been multiple results from combining Monte Carlo 
and DFT for such cases [34-37]. But none of these could 
approach the accuracy needed to invert the Kohn-Sham 
equations or extract highly accurate correlation energy 
components. For such purposes, finite molecular systems 
are often the only ones where sufficient accuracy can be 
practically achieved. 

The prototype case for electronic structure and chem¬ 
ical binding is the simplest molecule, H 2 , and its binding 
energy curve at zero temperature is simple to calculate, 
to study the success of GGA’s near equilibrium [38] and 
their failures as the bond is stretched [39]. But even this 
system is too difficult to calculate when the electrons are 
heated: Only the mean number of electrons is fixed, and 
all possible electron numbers must be included in evalu¬ 
ating the grand canonical partition function. 

Here we circumvent this difficulty with the simplest 
representation of a diatomic molecule. In a minimal ba¬ 
sis set (one function per atom), the full Hamiltonian is 
simply a 2-site Hubbard model to which lattice DFT 
applies [40]. The severe truncation of the Hilbert space 
makes exact solution possible in thermal DFT. By invert¬ 
ing the MKS equations, we perform the first exact cal- 
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culations of correlation free energies and their individual 
components for an inhomogeneous system, an admittedly 
crude representation of a chemical bond. By performing 
self-consistent calculations with the exact ground-state 
exchange-correlation energy functional for this system, 
we show that the ground-state approximation works well, 
even becoming relatively exact in the high-temperature 
limit. We also illustrate several exact conditions on the 
correlation energy components. While such a simplified 
model cannot be used to test the accuracy of standard 
approximations applied in the continuum, such as the 
local density or generalized gradient approximations, it 
does provide a first glimpse at the behaviors of correla¬ 
tion energy components as a function of temperature, a 
subject about which almost nothing is known outside of 
the uniform electron gas. 

This paper is laid out as follows. In section IIA we re¬ 
view the ground-state of the asymmetric Hubbard model. 
In section IIB we briefly outline thermal DFT. Next, in 
section III we write out the analytic expressions for the 
MB and MKS system. Lastly in section V we discuss 
some results using the ground-state XC functional. 

II. BACKGROUND 
A. Ground-state Hubbard Dimer 


very well characterized in the symmetric case, where an 
expansion in powers of U converges absolutely up to 
U = At and diverges beyond that; similarly, an expan¬ 
sion in 1/U converges absolutely only for U > At. Here, 
we restrict our attention to the weakly correlated regime 
in order to best mimic typical conditions of thermal DFT 
calculations. 


B. Thermal Density Functional Theory 

In this section we will briefly review the basics of ther¬ 
mal DFT [16]. For a more exhaustive treatment see Ref. 
[19]. We begin with an ensemble in thermal equilibrium 
connected to a bath at temperature r. The free energy 
may be found from: 

A = min ^F[n] + J d^rn(r)i;(r)^ (2) 

where i;(r) is the one-body potential, /i is the chemical 
potential, and the minimization is over all positive den¬ 
sities with finite kinetic energy. The Mermin functional 
is 

F[n] = mn Tr {(f + t>ee - T5)r} (3) 


Ref. [40] is an exhaustive review of the asymmetric 
Hubbard dimer for the ground-state case. In this section 
we briefly review the Hamiltonian and the most salient 
points. The Hamiltonian is typically written as 

H = (c\^C 2 o- + h.c)j H- y; {Uhi^hii + Vihi) (1) 

cr i 

where ct^(Qcr) are electron creation (annihilation) op¬ 
erators and Uia- = ct^Qcr are number operators, t is 
the strength of electron hopping between sites, U the 
Coulomb repulsion when two electrons are on the same 
site, and Vi is the external potential on each site. Without 
loss of generality, we choose vi^V2 = 0, /S.v = V2—Vi^ and 
denote the occupation difference An = n 2 —ni. All terms 
in Eq. (1) have analogs in an ab initio Hamiltonian [40]. 
The hopping term plays a role logically analogous to the 
kinetic energy, the Coulomb repulsion is now ultra-short 
ranged but otherwise the same, and the on-site potential 
serves as the one-body potential. Most importantly the 
asymmetry is necessary to perform our analysis. Other¬ 
wise the occupation difference would vanish and we could 
not learn about the function(al) behavior. We choose 
units where 2t = 1 and we vary U and /S.v. 

The key observation is that repulsion and asymmetry 
directly compete. When U dominates over /S.v the den¬ 
sity, An, tends towards 0, while in the opposite limit An 
tends towards 2. Additionally U < Av is the weakly- 
correlated regime while U > Av is strongly-correlated. 
The difference between weak and strong correlation is 


where T is the kinetic energy^ operator, We the electron- 
electron repulsion operator, S the entropy operator, and 
the minimization is over all statistical density matrices 
with density n(r). The average particle number is deter¬ 
mined uniquely by ji. Then one can construct the MKS 
equations [19] 

+ <[n](r)| 4>1{y) = ej 4>1 {y), (4) 


where 


< [n] (r) = (r) + Vh [n] (r ) + [n] (r) , (5) 

and Vu[n]{r) is just the usual Hartree potential[41] and 


i„(r) ■ 

The density is the sum over all orbitals. 


(6) 


( 7 ) 


where fj = are their Fermi occupations. 

Finally, once self-consistency has been achieved, the free 
energy of the interacting system is reconstructed as: 

= As - C/h [n] + AJc [n] - J M (’’) • (8) 

where Ag is the Kohn-Sham free energy. 



3 


If the exact XC free energy density functional (confus¬ 
ingly, often referred to as simply the XC energy) were 
known and used in the MKS equations, then their solu¬ 
tion produces the exact density and free energy (and any 
other quantity that can be directly extracted from them). 
However, there are very few cases where we have access 
to the exact 'r’xc(i*)- All practical MKS calculations use 
some approximation, and most use a simple ground-state 
approximation. To distinguish different levels of approx¬ 
imation, we write 

A-xcW ~ (9) 

where ^xc['^] is the exact ground-state XC energy, and 
AA^^[n] is the difference in XC free energy from its 
ground-state value. We call this the thermal contribu¬ 
tion to Then, the zero-temperature approximation 

(ZTA) is where we ignore the thermal contribution to 
4^ i P 

^XC5 


where Ei{N) is the Tth energy level of the Hamiltonian 
with N particles. The energies for 0 through 4 parti¬ 
cles are calculated explicitly, yielding the exact partition 
function. From that we construct the grand potential, its 
derivatives, and the free energy in the usual fashion[43]: 


n = -T\Og{Zgc), 




A = fiN — Q. 


( 12 ) 

(13) 


We choose half-filling, {N) = 2, which means ja = U{2 
(and /i = 0 for the MKS system) [44]. With the partition 
function and Boltzmann factors we can calculate ensem¬ 
ble averages:. 


X = (14) 

i,N 


= E^c[n\, ( 10 ) 

i.e., we neglect thermal effects, but use the exact ground- 
state XC functional. This allows us to separate thermal 
from non-thermal XC effects in a completely well-defined 
manner. Of course, in practice, it is only in simple model 
systems that one has access to the exact ground-state XC 
functional. 

In this language, most modern QMD calculations can 
be thought to have made two distinct approximations. 
The first is to make ZTA and ignore thermal contribu¬ 
tions. The second is to use some common approximation 
for F^xc['^] within ZTA. On the other hand, calculations 
that use, e.g., thermal LDA, go beyond ZTA, but approx¬ 
imate both the ground-state and thermal contributions 
to Alg.[n\. 

III. ANALYTIC RESULTS 

We apply this technology to the asymmetric Hubbard 
dimer. The DFT version of a lattice model is called site- 
occupation functional theory (SOFT) [42] and has the dis¬ 
tinct advantage of a truncated Hilbert space. We can 
compute every energy for every particle number and con¬ 
struct exact thermodynamic and DFT components. The 
truncation makes the calculation feasible. We expect 
that, for very high temperatures, the results will not be 
representative of realistic systems with infinite Hilbert 
spaces. 

A. Exact many-body solution 

To begin, we calculate the finite-temperature many- 
body energy and density for the Hubbard dimer. Begin 
with the grand canonical partition function 

= ( 11 ) 

i,N 


where {X)i^N is the expectation value of a general oper¬ 
ator X of the i-th state for N particles. Using Eq. (14) 
we compute the exact energy components for the dimer. 
To do this, we calculate the expectation values for each 
particle number of the quantities of interest such as T, 
Uee, and An. We list in the appendix all the expectation 
values for the total energies, energy components, coeffi¬ 
cients of the eigenstates, and densities for all the particle 
numbers. 



T 


FIG. 1. Free energy for different values of Av. Solid lines are 
exact, dashed lines are the zero-temperature XC approxima¬ 
tion (ZTA), evaluated on the self-consistent thermal density. 

In Fig. 1 and 2, we plot the free energy and entropy as 
a function of temperature for several different values of 
Av. For the free energy we include curves for the zero- 
temperature approximation and for the entropy we in¬ 
clude the self-consistent Kohn-Sham entropy (both to be 
discussed later). In both cases we pick a system, i.e. fix 
Av and U and see what happens as we heat it up. For the 
free energy, the values at r = 0 recover the ground-state 
energies reported in Ref. [40]. Increasing temperature 
results in a decrease in free energy primarily due to the 
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FIG. 2. Exact entropy (solid) and self-consistent Kohn-Sham 
entropy (dashed) for different values of Av. All curves ap¬ 
proach 4 log 2. 

entropic term, —r5', as expected. At small temperatures 
there is minimal effect as seen in Fig. 2 where the entropy 
is small and further multiplied by a r ^ 1 when calculat¬ 
ing A. However, once the system is sufficiently warm the 
entropy plays a much larger role. In contrast, increas¬ 
ing Av lowers the entropy since the asymmetry restricts 
the motion of electrons. Lastly, the entropy approaches a 
maximum value of log(16) for higher temperatures where 
16 is the number of states in our grand canonical ensem¬ 
ble. 


B. Inversion and correlation components 

Next, we construct the exact KS potential as well as 
various energy components using the MKS approach. To 
begin we construct the exact occupation difference An 
from Eq. (14). We plot the result in Fig. 3 for fixed U 
but against Av and vary r. In this figure we also plot 
the ZTA result which will be discussed later. Increasing 
the temperature pushes the electrons apart, akin to re¬ 
pulsion. As the system heats up. An becomes closer to 
0 as both electrons sit on separate sites even when Av is 
large. 

To construct the exact MKS potential, we first give for¬ 
mulas for non-interacting electrons {U = 0, a.k.a. tight- 
binding). 

The grand canonical partition function collapses to the 
product 

= + (15) 

where is the single-particle orbital energy. Eq. (12) 
and (13) can then be used. The entropy can also be 
explicitly given in terms of Fermi factors, 

Ss = -'^fi log(/i) + (1 - fi) log(l - fi). (16) 



FIG. 3. Densities as a function of temperature for the system 
of Fig. 1. Solid lines are exact, dashed lines are self-consistent 
KS using the ZTA. 

where fi = (1 + Kohn-Sham entropy 

is calculated in Fig. 2 where the Fermi factors are cal¬ 
culated from self-consistently solving the MKS equations 
(see below). 

To construct the MKS system for the Hubbard dimer 
within SOFT, we simply repeat the exact calculation 
with [/ = 0, i.e., a tight-binding dimer. We find: 

An = — 2 sin^tanha (17) 

where a = (4r cos0)“^, smcj) = x/\J\ + r is in units 
of 2t, and x = Aug/2t. To perform the inversion for a 
given density from the many-body problem, we perform 
a binary search at the given temperature on Eq. (17) 
to find Aus(An), the exact KS site-potential difference 
that yields the required occupation density. The exact 
Auxc for fho given An is then found by subtracting off 
the other potential contributions, i.e., Av and Auh- The 
Hartree energy (in the standard DFT definition [41]) for 
this model is 

f/H(An) = u(l+ , (18) 

and the Hartree potential is simply 

UH(An) = U An/2 (19) 

and both functionals are temperature-independent. For 
two unpolarized electrons, Ey^ = —Uii/2 at all 
temperatures [19], and so is also independent of r. The 
thermal MKS hopping energy is just that of this tight- 
binding problem: 

Tf {An)/{2t) = An/x{An) (20) 

and the tight-binding MKS entropy is 

Sf {An) = 4 log {2 cosh a} — Aa tanh a (21) 

With these simple results, we can now extract the corre¬ 
lation free energy for this problem as 

= (T- - n) - riS^ - SI) + (v; - C^Hx) (22) 
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where 5"^, and are calculated from the many- 
body problem via eqs. (14), (12), and (14). Since is 
trivial and has no thermal contribution for our system, 
is what we study, and we know of no other exact 
calculation of this quantity for a finite system. 


IV. NUMERICAL RESULTS 

Performing the inversion to explicitly analyze the MKS 
potential shows how the features of interactions are built 
into the non-interacting potential [45-47]. The crux of 
the MKS approach is that we capture the effects of in¬ 
teractions through the modified external potential Avq. 
For example, interaction causes the dimer occupations 
to be more symmetric, thus Avq < Av for a MB sys¬ 
tem with U > 0. Similarly, for any given density both 
potentials, Av and Avs, increase with temperature to 
counteract thermal effects pushing the system towards 
symmetry. But even in this simple model, there is a vast 
parameter space to be explored as, choosing 2t = 1, we 
can vary U, Av, r, and (N). We focus on {N) = 2, 
and the weakly-correlated and low temperature corner of 
our parameter space: U,r < 1. In particular, we avoid 
warming our model so much that properties are strongly 
influenced by the very limited Hilbert space. Specifically, 
we check that the system is not too hot by computing the 
occupations of all the states in the grand canonical en¬ 
semble. We test this in the symmetric case because it 
is most prone to overheating since asymmetry competes 
against thermal effects. For U = 1, uniform occupation 
of all states does not occur until r S and appreciable 
uniformity does not start to arise until r ^ 4. Thus our 
results are not limited by the top of our Hilbert space. 

We can calculate all the individual contributions to the 
correlation free energy by subtracting MKS quantities 
from their physical counterparts. These are the energy 
differences appearing in Eq. (22): 

Tr=T-_T-^ S- = S^-Sl, UJ = K;-Uhx. (23) 

The kentropic correlation is — r Sc and plays a 

key role in thermal DFT[48]. In Fig. 4, we plot the exact 
correlation free energy functional, the sum of kinetic and 
potential correlation functional, and lastly the entropic 
correlation functional all for various temperatures. By 
fixing U and r and plotting versus An, we analyze the 
correlation as a density functional, i.e. we are no longer 
looking at a fixed system and instead are looking at the 
underlying structure of how thermal DFT behaves. 

We see that the correlation free energy is always 
negative, the kentropic contribution is always positive 
(not shown), and the potential contribution is always 
negative. These are consistent with conditions on the 
correlation [48]. This is the first exact investigation of 
those inequalities. The correlation free energy, al¬ 
ways decreases with temperature at U = 1, even though 
the components do not behave that way at small tem¬ 
perature. Tc + Uc and rSc also decrease for all densities 





An 


FIG. 4. Panel 1: Correlation free energy functional for 
various temperatures. Panel 2: Sum of kinetic and potential 
energy functional for various temperatures. Panel 3: Entropic 
correlation functional for various temperatures. 


at larger temperature just like Ac. In this regime, ther¬ 
mal effects dominate over interactions, resulting in the 
interacting system and the non-interacting system hav¬ 
ing similar energy components and thus relatively smaller 
correlation. But for small temperature, i.e. r < 1 when 
U = 1, the MKS quantities are furthest from the exact 
system since neither effect dominates and this results in 
an even larger difference between the two systems than 
at r = 0. Overall we see the same behavior as in the 
ground-state case [40] - correlation decreases as our sys¬ 
tem becomes more asymmetric. If the electrons are com¬ 
pletely pinned on the lower site then there is no motion. 
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the interaction is completely described by the Hartree, 
and there is only one entropic conformation. 



A 

FIG. 5. Adiabatic connection integrand for the symmetric 
dimer at several different temperatures. 



In Fig. 5, for the symmetric case, turning on temperature 
clearly reduces both the magnitude of the correlation and 
the degree of static correlation, as judged by the initial 
slope of the curves. Fig. 6 shows this result remains true 
beyond the symmetric case. 



T 

FIG. 7. Gorrelation free energy for the symmetric case with 
increasing values of U ranging from weak to strong correla¬ 
tion. 

In Fig. 7, we repeat the curves of Fig. 4 but now 
for fixed An = 0 and increasing U. We start with the 
U = 1 from earlier and increase into the strongly cor¬ 
related regime. The curves show a minimum at about 
r = 0.25, particularly in U = 3 and 4. Thus the deriva¬ 
tive with respect to temperature can be negative, and this 
does not happen even if we look closely nt U = 1. Thus 
the correlation free energy is not generally monotonically 
decreasing in magnitude and the correlation energy is not 
bounded by the r = 0 value. 


V. ZERO-TEMPERATURE APPROXIMATION 


FIG. 6. Same as Fig. 5, except An — 1. 

Next, we consider the adiabatic connection formula [49, 
50] that has proven useful in studying and improving den¬ 
sity functional approximations. The ground-state ver¬ 
sion was calculated for the Hubbard dimer in Fig. 21 
of Ref. [40]. An alternative version, called the thermal 
connection formula, was derived in Ref. [51], but that fia- 
vor relies on relating the coupling-constant to coordinate 
scaling. Such a procedure applies to continuum models, 
but not lattices. So we use the traditional version here, 
applied to finite temperature [48]: 



where A is a coupling constant inserted in front of We in 
the Hamiltonian, but (unlike regular many-body theory) 
the density is held fixed during the variation. Here Uc^ 
is the potential correlation energy at coupling constant A, 
which, for our model, is obtained by replacing U with XU. 


In this section, we explore the effects of making the 
zero-temperature approximation (ZTA), in which ther¬ 
mal contributions are ignored (Eq. (10)). We use the 
(essentially) exact parametrization of the ground-state 
XC energy of the Hubbard dimer of Eq. (108) of Ref. 
[40]. This substitution is made in the calculation of the 
total free energy and in the MKS equations via the cal¬ 
culation of the XC potential, Eq. (6). We return to 
Eig. 1, where we also plot the free energy in the ZTA by 
replacing A^(An) with Ec(An), evaluated on the self- 
consistent An. We see that the error of ZTA is ex¬ 
tremely small for r < 0.5. Moreover, trends are very 
well reproduced by the ZTA values, and fractional er¬ 
rors shrink for large r. This suggests that free energies 
in such calculations may be reliable depending, of course, 
on the precision needed in a given calculation. The errors 
grow most rapidly with r when the dimer is asymmet¬ 
ric. Thermal effects push the electrons apart, making the 
density more symmetric, in direct competition with Av. 
Eor larger An, there is a larger error in ignoring ther¬ 
mal effects. Note that since we have only two electrons. 
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our model is a worst case scenario. In many simulations, 
there are more valence electrons per site, and (exchange- 
)correlation components are a much smaller fraction of 
the total energy. In a realistic DFT calculation, the er¬ 
ror made by approximating the ground-state functional 
would likely be much larger than the error due to the 
lack of temperature-dependence [52]. 



0 1 2 3 4 5 


Av 

FIG. 8. Error in ZTA densities of Fig. 3, density from self- 
consistent MKS subtracted from exact density. 

However, this is only part of the story. Real thermal 
DFT calculations are performed self-consistently within 
ZTA. Then both the density and MKS orbitals are of¬ 
ten used to calculate response properties (usually on the 
MKS orbitals) [53-61]. In Fig. 3, we compare the self- 
consistent density obtained using Eq. (4) through Eq. 
(7). In Fig. 8, we plot the differences. We see that the 
maximum errors in the density are small. At first they 
grow with small temperature but quickly start to lessen 
as temperature increases which will be further explained 
below. As Av gets large the error goes to zero since the 
asymmetry dominates over thermal effects. 

In terms of Fig. 4, the ZTA consists of approximat¬ 
ing each of the curves by the corresponding black one. 
Because all correlation components tend to vanish with 
increasing temperature, while the total free energy grows 
in magnitude, the small error made in the ZTA becomes 
less relevant with increasing temperature. Specifically, 
we can analyze the symmetric case where correlation ef¬ 
fects are at their strongest. At r = 0 correlation is about 
20% of the total energy but when the system is at r = 1 
correlation is roughly 2.5% of the total free energy. More 
importantly this is due to the total energy magnitude 
going up by a factor of 5 and the correlation only de¬ 
creasing by a factor of 2. This explains the small errors 
in the ZTA free energies of Fig. 1 and the behavior of 
the self-consistent ZTA densities of Fig. 8. Note that the 
temperatures need not be so high as to make the density 
uniform (i.e., symmetric). Fig. 3 shows that, even for the 
temperature at which density differences can be largest 
(r = 1), the density difference can remain substantial as 
the temperature increases, if the inhomogeneity (An) is 
large enough. 


VI. CONCLUSIONS 

In summary, we have solved the simplest possible non¬ 
trivial system at finite temperature exactly, both for the 
many-body case and within MKS density functional the¬ 
ory. We have produced the first exact plots of MKS 
quantities and the ZTA approximation for a finite sys¬ 
tem (albeit one with a limited Hilbert space). When the 
system is weakly correlated system at low to moderate 
temperatures, the neglect of thermal contributions to the 
exchange-correlation functional has relatively little effect 
on the calculated free energies and even less on the self- 
consistent densities. Present limitations of ground-state 
approximations, such as their inability to treat strongly 
correlated systems, are likely the greatest source of er¬ 
ror in these calculations. Future work will explore other 
quantities of interest within thermal DFT and will ana¬ 
lyze the ZTA more deeply. 
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Appendix: Energies and Densities for all States 

Here we list all the total energies, energy components, 
and density components for all particle numbers so that 
all the relevant ensemble averages of Eq. (14) can be 
reconstructed. We begin with the energies 

Ei{A -N) = {2- N)U + Ei{N) N = 0-4, 

^o(O) = 0, 

Eo,iV) = TV(2i)2 + Av2/2, 

2TI 2r 

^i(2) = -^ cos(6> + —(i + 1)) i = 0,1, 2, 

F;,(2) =0, 2 = 3,4, 5, 

where 

r = y/3((2t)^ -|- Au^) + 

^ 1 r 9U(Av‘^ -U^ 

0-- arccos [( 3 (( 2 t )2 + Ai; 2 ) + C/2)3/2_ ' 

Ei{2) and E 2 { 2 ) are both positive and should be ordered 
4 and 5 instead. However the three triplets, i.e. the three 












zero-energy states, give only zero values in the later ex¬ 
pectation values, so for notational convenience we order 
the non-zero 2-particle states 0, 1, and 2 instead of 0, 4, 
and 5. These energies were used to construct 'Z^gc in Eq. 

(11) of the main text. 

Next are the expectation values needed to construct 
the three different ensemble averages of interest, T, Eee, 
and An (Eext is unnecessary since it is trivially Ai;An/2): 


Ti(4-A) =Ti(n) A = 0-4 

r,(2) = (/3++/3-)VE,(2) i = 0,l,2, 

Ve,i(4 - A^) = (2 - A^)[/ + V;e,i(iV) N = 0-4, 
Ve,0,l(l) =0, 

Ve,i(2) = U{{^tf + (/?*■)') ^ = 0,1,2, 


Ani(4-Ar) = Ani(A^) A = 0-4, 

A /i\ 2 Av 

Ano i(l) = =F , , 

y(2i)2 + Av2’ 

Ani(2)=2((^r)2_(^+)2) 


i = 0,1,2, 


and all the 0-particle terms are 0. The ;d^’s are from the 
A = 2 wavefunction: 


\MN)) = ai(iV)(|12) + |21)) +/3+(A)|ll) + ^r(A)|22) 


with 

^2t{Ei{2)-U) ^^U-Ei{2)±Av 

CiEi{2) ’ Ci 

Ci = v/2(At;2 + {Ei{2) - U^il + {2t/Ei{2)Y). 


The ket \ij) signifies an electron at site i and site j. These 
expectation values were used with Eq. (14) to construct 
the densities and energy components shown in the figures. 


[1] P. Hohenberg and W. Kohn, “Inhomogeneous electron 
gas,” Phys. Rev. 136, B864-B871 (1964). 

[2] Frank Graziani, Michael P. Desjarlais, Ronald Redmer, 
and Samuel B. Trickey, eds., Frontiers and Challenges in 
Warm Dense Matter^ Lecture Notes in Computational 
Science and Engineering, Vol. 96 (Springer International 
Publishing, 2014). 

[3] U.S. Department of Energy, Basic Research Needs for 
High Energy Density Laboratory Physics: Report of the 
Workshop on High Energy Density Laboratory Physics 
Research Needs, Tech. Rep. (Office of Science and Na¬ 
tional Nuclear Security Administration, 2009). 

[4] Thomas R. Mattsson and Michael P. Desjarlais, “Phase 
diagram and electrical conductivity of high energy- 
density water from density functional theory,” Phys. Rev. 
Lett. 97, 017801 (2006). 

[5] Winfried Lorenzen, Bastian Holst, and Ronald Redmer, 
“Demixing of hydrogen and helium at megabar pres¬ 
sures,” Phys. Rev. Lett. 102, 115701 (2009). 

[6] M. D. Knudson, M. P. Desjarlais, A. Becker, R. W. 
Lemke, K. R. Cochrane, M. E. Savage, D. E. Bliss, T. R. 
Mattsson, and R. Redmer, “Direct observation of an 
abrupt insulator-to-metal transition in dense liquid deu¬ 
terium,” Science 348, 1455-1460 (2015). 

[7] M. D. Knudson and M. P. Desjarlais, “Shock compression 
of quartz to 1.6 TPa: Redefining a pressure standard,” 
Phys. Rev. Lett. 103, 225501 (2009). 

[8] M. D. Knudson, M. P. Desjarlais, and A. Pribram-Jones, 
“Adiabatic release measurements in aluminum between 
400 and 1200 gpa: Characterization of aluminum as a 
shock standard in the multimegabar regime,” Phys. Rev. 
B 91, 224105 (2015). 

[9] Bastian Holst, Ronald Redmer, and Michael P. Desjar¬ 
lais, “Thermophysical properties of warm dense hydrogen 
using quantum molecular dynamics simulations,” Phys. 
Rev. B 77, 184201 (2008). 

[10] Andre Kietzmann, Ronald Redmer, Michael P. Desjar¬ 
lais, and Thomas R. Mattsson, “Complex behavior of 
fluid lithium under extreme conditions,” Phys. Rev. Lett. 


101, 070401 (2008). 

[11] Seth Root, Rudolph J. Magyar, John H. Carpenter, 
David L. Hanson, and Thomas R. Mattsson, “Shock 
compression of a fifth period element: Liquid xenon to 
840 GPa,” Phys. Rev. Lett. 105, 085501 (2010). 

[12] R F Smith, J H Eggert, R Jeanloz, T S Duffy, D G Braun, 
J R Patterson, R E Rudd, J Biener, A E Lazicki, A V 
Hamza, J Wang, T Braun, L X Benedict, P M Celliers, 
and G W Collins, “Ramp compression of diamond to five 
terapascals,” Nature 511, 330-3 (2014). 

[13] MWC Dharma-Wardana and R Taylor, “Exchange and 
correlation potentials for finite temperature quantum 
calculations at intermediate degeneracies,” Journal of 
Physics C: Solid State Physics 14, 629 (1981). 

[14] Francois Perrot and M. W. C. Dharma-wardana, “Ex¬ 
change and correlation potentials for electron-ion sys¬ 
tems at finite temperatures,” Phys. Rev. A 30, 2619-2626 
(1984). 

[15] John P. Perdew, Kieron Burke, and Matthias Ernzer- 
hof, “Generalized gradient approximation made simple,” 
Phys. Rev. Lett. 77, 3865-3868 (1996), ibid. 78, 1396(E) 
(1997). 

[16] N. D. Mermin, “Thermal properties of the inhomogenous 
electron gas,” Phys. Rev. 137, A: 1441 (1965). 

[17] W. Kohn and L. J. Sham, “Self-consistent equations in¬ 
cluding exchange and correlation effects,” Phys. Rev. 
140, A1133-A1138 (1965). 

[18] S. Atzeni and J. Meyer-ter Vehn, The Physics of Inertial 
Fusion: Beam-Plasma Interaction, Hydrodynamics, Hot 
Dense Matter (Clarendon Press, 2004). 

[19] Aurora Pribram-Jones, Stefano Pittalis, E.K.U. Gross, 
and Kieron Burke, “Thermal density functional theory 
in context,” in Frontiers and Challenges in Warm Dense 
Matter, Lecture Notes in Computational Science and En¬ 
gineering, Vol. 96, edited by Frank Graziani, Michael P. 
Desjarlais, Ronald Redmer, and Samuel B. Trickey 
(Springer International Publishing, 2014) pp. 25-60. 

[20] Valentin V. Karasiev, Travis Sjostrom, and S. B. 
Trickey, “Comparison of density functional approxima- 










9 


tions and the finite-temperature hartree-fock approxima¬ 
tion in warm dense lithium,” Phys. Rev. E 86, 056704 
( 2012 ). 

[ 21 ] M. Ernzerhof and G. E. Scuseria, “Assessment of 
the Perdew-Burke-Ernzerhof exchange-correlation func¬ 
tional,” J. Chem. Phys 110, 5029 (1999). 

[22] M. Marsman J. Paier, R. Hirschl and G. Kresse, “The 
perdew-burke-ernzerhof exchange-correlation functional 
applied to the g 2 -l test set using a plane-wave basis set,” 
J. Ghem. Phys. 122, 234102 (2005). 

[23] P.R.T. Schipper O.V. Gritsenko and E.J. Baerends, “Ex¬ 
change and correlation energy in density functional the¬ 
ory. comparison of accurate dft quantities with tradi¬ 
tional hartree-fock based ones and generalized gradient 
approximations for the molecules 62 , n 2 , and f 2 ,” J. 
Ghem. Phys. 107, 5007 (1997). 

[24] Paul Ziesche, Stefan Kurth, and John P. Perdew, “Den¬ 
sity functionals from LDA to GGA,” Gomputational Ma¬ 
terials Science 11, 122 - 127 (1998). 

[25] J.P. Perdew S. Kurth and P. Blaha, “Molecular and solid- 
state tests of density functional approximations: Lsd, 
gga’s, and meta-gga’s,” Int. J. Quantum Ghem. 75, 889 
(1999). 

[26] Anubhav Jain, Shyue Ping Ong, Geoffroy Haulier, 
Wei Ghen, William Davidson Richards, Stephen Dacek, 
Shreyas Gholia, Dan Gunter, David Skinner, Gerbrand 
Geder, and Kristin A. Persson, “Gommentary: The ma¬ 
terials project: A materials genome approach to acceler¬ 
ating materials innovation,” APL Materials 1, - (2013). 

[27] V S Eilinov, M Bonitz, W Ebeling, and V E Eortov, 
“Thermodynamics of hot dense h-plasmas: path integral 
monte carlo simulations and analytical approximations,” 
Plasma Physics and Gontrolled Eusion 43, 743 (2001). 

[28] T. School, M. Bonitz, A. Eilinov, D. Hochstuhl, and 
J.W. Dufty, “Gonfiguration path integral monte carlo,” 
Gontributions to Plasma Physics 51, 687-697 (2011). 

[29] K. P. Driver and B. Militzer, “All-electron path integral 
monte carlo simulations of warm dense matter: Appli¬ 
cation to water and carbon plasmas,” Phys. Rev. Lett. 
108, 115502 (2012). 

[30] Ethan W. Brown, Bryan K. Glark, Jonathan L. DuBois, 
and David M. Geperley, “Path-integral monte carlo sim¬ 
ulation of the warm dense homogeneous electron gas,” 
Phys. Rev. Lett. 110, 146405 (2013). 

[31] T. School, S. Groth, J. Vorberger, and M. Bonitz, “A 6 
Initio thermodynamic results for the degenerate electron 
gas at hnite temperature,” Phys. Rev. Lett. 115, 130402 
(2015). 

[32] Burkhard Militzer and Kevin P. Driver, “Development 
of path integral monte carlo simulations with localized 
nodal surfaces for second-row elements,” Phys. Rev. Lett. 
115, 176403 (2015). 

[33] S. Groth, T. School, T. Dornheim, and M. Bonitz, “A 6 
initio quantum monte carlo simulations of the uniform 
electron gas without fixed nodes,” Phys. Rev. B 93, 
085102 (2016). 

[34] B. Militzer, “Path integral monte carlo and density func¬ 
tional molecular dynamics simulations of hot, dense he¬ 
lium,” Phys. Rev. B 79, 155105 (2009). 

[35] J. Vorberger and D. O. Gericke, “A 6 initio approach to 
model x-ray diffraction in warm dense matter,” Phys. 
Rev. E 91, 033112 (2015). 

[36] K. P. Driver and B. Militzer, “Eirst-principles simulations 
and shock hugoniot calculations of warm dense neon,” 


Phys. Rev. B 91, 045103 (2015). 

[37] K. P. Driver and B. Militzer, “Eirst-principles equation 
of state calculations of warm dense nitrogen,” Phys. Rev. 
B 93, 064101 (2016). 

[38] David G. Patton, Dirk V. Porezag, and Mark R. Peder¬ 
son, “Simplified generalized-gradient approximation and 
anharmonicity: Benchmark calculations on molecules,” 
Phys. Rev. B 55, 7454-7459 (1997). 

[39] E.J. Baerends, “Exact exchange-correlation treatment of 
dissociated h .2 in density functional theory,” Phys. Rev. 
Lett. 87, 133004 (2001). 

[40] D J Garrascal, J Eerrer, J G Smith, and K Burke, 
“The hubbard dimer: a density functional case study of 
a many-body problem,” Journal of Physics: Gondensed 
Matter 27, 393001 (2015). 

[41] A. Pribram-Jones, Z.-H. Yang, J. R. Trail, K. Burke, 
R. J. Needs, and G. A. Ullrich, “Excitations and bench¬ 
mark ensemble density functional theory for two elec¬ 
trons,” J. Ghem. Phys. 140, 18A541 (2014). 

[42] K. Schonhammer, O. Gunnarsson, and R. M. Noack, 
“Density-functional theory on a lattice: Gomparison with 
exact numerical results for a model with strongly corre¬ 
lated electrons,” Phys. Rev. B 52, 2504-2510 (1995). 

[43] Paul D. Beale, R.K. Pathria, Statistical Mechanics, 3rd 
ed. (Academic Press, 2011). 

[44] Hiroyuki Shiba, “Thermodynamic properties of the one¬ 
dimensional half-filled-band hubbard model, ii applica¬ 
tion of the grand canonical method,” Progress of Theo¬ 
retical Physics 48, 2171-2186 (1972). 

[45] M. Thiele, E. K. U. Gross, and S. Kiimmel, “Adia¬ 
batic approximation in nonperturbative time-dependent 
density-functional theory,” Phys. Rev. Lett. 100, 153004 
(2008). 

[46] J. P. Goe, K. Gapelle, and 1. D’Amico, “Reverse engi¬ 
neering in many-body quantum physics: Gorrespondence 
between many-body systems and effective single-particle 
equations,” Phys. Rev. A 79, 032504 (2009). 

[47] S. E. B. Nielsen, M. Ruggenthaler, and R. van Leeuwen, 
“Many-body quantum dynamics from the density,” EPL 
(Europhysics Letters) 101, 33001 (2013). 

[48] S. Pittalis, G. R. Proetto, A. Eloris, A. Sanna, G. Bersier, 
K. Burke, and E. K. U. Gross, “Exact conditions in 
hnite-temperature density-functional theory,” Phys. Rev. 
Lett. 107, 163001 (2011). 

[49] D.G. Langreth and J.P. Perdew, “The exchange- 
correlation energy of a metallic surface,” Solid State 
Gommun. 17, 1425 (1975). 

[50] O. Gunnarsson and B.I. Lundqvist, “Exchange and corre¬ 
lation in atoms, molecules, and solids by the spin-density- 
functional formalism,” Phys. Rev. B 13, 4274 (1976). 

[51] “Gonnection formula for thermal density functional the¬ 
ory,” Phys. Rev. (2015), submitted. 

[52] Travis Sjostrom and Jerome Daligault, “Gradient correc¬ 
tions to the exchange-correlation free energy,” Phys. Rev. 
B 90, 155109 (2014). 

[53] Pier Luigi Silvestrelli, “No evidence of a metal-insulator 
transition in dense hot aluminum: A first-principles 
study,” Phys. Rev. B 60, 16382-16388 (1999). 

[54] Pier Luigi Silvestrelli, Ali Alavi, and Michele Parrinello, 
“Electrical-conductivity calculation in ab initio simula¬ 
tions of metals:application to liquid sodium,” Phys. Rev. 
B 55, 15515-15522 (1997). 

[55] Vanina Recoules, Patrick Renaudin, Jean Glerouin, 
Pierre Noiret, and Gilles Zerah, “Electrical conductivity 


10 


of hot expanded aluminum: Experimental measurements 
and ab initio calculations,” Phys. Rev. E 66, 056412 
( 2002 ). 

[56] Jean Clerouin, Patrick Renaudin, Yann Laudernet, 
Pierre Noiret, and Michael P. Desjarlais, “Electrical con¬ 
ductivity and equation-of-state study of warm dense cop¬ 
per: Measurements and quantum molecular dynamics 
calculations,” Phys. Rev. B 71, 064203 (2005). 

[57] M. P. Desjarlais, J. D. Kress, and L. A. Collins, “Electri¬ 
cal conductivity for warm, dense aluminum plasmas and 
liquids,” Phys. Rev. E 66, 025401 (2002). 

[58] Vanina Recoules, Jean Clerouin, P Renaudin, P Noiret, 
and Gilles Zrah, “Electrical conductivity of a strongly 
correlated aluminium plasma,” Journal of Physics A: 


Mathematical and General 36, 6033 (2003). 

[59] M. P. Desjarlais, “Density functional calculations of the 
reflectivity of shocked xenon with ionization based gap 
corrections,” Contributions to Plasma Physics 45, 300- 
304 (2005). 

[60] Vanina Recoules, Elavien Lambert, Alain Decoster, 
Benoit Canaud, and Jean Clerouin, “A6 Initio determi¬ 
nation of thermal conductivity of dense hydrogen plas¬ 
mas,” Phys. Rev. Lett. 102, 075002 (2009). 

[61] K.-U. Plagemann, P. Sperling, R. Thiele, M. P. Desjar¬ 
lais, C. Eortmann, T. Doppner, H. J. Lee, S. H. Glen- 
zer, and R. Redmer, “Dynamic structure factor in warm 
dense beryllium,” New Journal of Physics 14, 055020 
( 2012 ). 


